clear
cd ..\..

local fname analyze_2_figures
capture log close
local c_date = c(current_date)
local c_time = c(current_time)
local c_time_date = "`c_date'"+"_"+"`c_time'"
local c_time_date = subinstr("`c_time_date '", ":", "_", .)
local c_time_date = subinstr("`c_time_date '", " ", "_", .)
log using "output/logs/`fname'_`c_time_date'.txt", replace text

********************************************************************************
*                            Define Macros                                     *
******************************************************************************** 
macro drop _all

* specify baseline consumption variable
global pre_avg "pre_avg_12"

* common covariates
global b $pre_avg // inherit baseline consumption variable from pre_avg
global c "income member building_size uni_size vintage married"
global waveXmonth "month#opower_paper_wave"

********************************************************************************
*                         Opower wave summary                                  *
******************************************************************************** 
use .\intermediate_data\elec_2014-2018_mainfile, clear

* collapse data, add wave 5 (gas only), and remove missing waves
bys ky_ba: gen opower_unique_sep = 1 if _n==1
collapse (sum) opower_unique_sep (mean) tot_kwh, by(opower_paper_wave opower)
tsset opower opower_paper_wave
tsfill
drop if opower_paper_wave==. | opower==.

graph bar (mean) opower_unique_sep, over (opower) over (opower_paper_wave) asyvars ///
	title("Number of Electric Accounts by Opower Wave") ///
	ytitle ("Number of Accounts") ylabel(, format(%9.0fc)) ///
	legend(label(1 "Control") label(2 "Treatment")) ///
	graphregion(color(white)) ///
	name(household_wave, replace)

replace tot_kwh = 0 if tot_kwh==. // wave 5
replace opower_unique_sep = 1 if opower_unique_sep==. // wave 5
graph bar (mean) tot_kwh [fweight=opower_unique_sep], over (opower_paper_wave) ///
	title("Average Consumption by Opower Wave") ///
	ytitle("Average Monthly Electricity Consumption (kwh)") ///
	graphregion(color(white)) ///
	name(mean_kwh_wave, replace)

graph combine household_wave mean_kwh_wave, ///
	cols (1) ysize(12) xsize(12) graphregion(color(white)) altshrink

graph export .\output\figures\FC.pdf, replace


********************************************************************************
*                             Data trend                                       *
******************************************************************************** 
use .\intermediate_data\elec_wave367_winterpre.dta, clear

foreach i of numlist 3 6 7 {
	binscatter tot_kwh month_since_opower if opower_paper_wave==`i' & month_since_opower<=24 & month_since_opower>=-24, ///
		discrete by(opower) line(connect) ///
		colors(gs4 blue) msymbol(D 0) ///
		xlabel(-24(12)24) ylabel(300(100)800, nogrid) ///
		xtitle("Months since Opower Assignment") ytitle("Monthly Electricity Consumption (kwh)") ///
		legend(label(1 "Control") label(2 "Treatment")) ///
		xline(0, lpattern(dash)) name(F2_wave`i', replace)
	graph export .\output\figures\F2_wave`i'.pdf, replace
}

********************************************************************************
*                            Event study graphs                                *
******************************************************************************** 

// wave3 graph------------------------------------------------------------------
use if opower_paper_wave==3 using .\intermediate_data\elec_wave367_winterpre.dta, clear

*bundle event month outside of [-22,21] range
replace month_since_opower=-23 if month_since_opower<-22
replace month_since_opower=22 if month_since_opower>21
* generate event time dummy
tab month_since_opower, gen (m_o)
* generate opower_event time interaction dummy
foreach i of numlist 1/25 {
	gen int_month_since_opower`i'=m_o`i'*opower 
}
* rename all event month before opower as month_o_i
foreach i of numlist 1/2 {
	local j=3-`i'
	rename int_month_since_opower`i' month_o_`j' 
}
* rename all event month after opower as month_oi
foreach i of numlist 4/25 {
	local j=`i'-3
	rename int_month_since_opower`i' month_o`j'
}
* rename start month of opower as month_o
rename int_month_since_opower3 month_o

// run regression
statsby _b _se, clear: reghdfe tot_kwh month_o_2 month_o month_o1-month_o22 $b $c, ///
	absorb($waveXmonth) vce(cluster ky_ba)

* generate coefficient and confidence intervals for each event month
gen o_coef_1=_b_month_o_2
gen o_se_1=_se_month_o_2
gen o_coef_2=0
gen o_se_2=0
gen o_coef_3=_b_month_o
gen o_se_3=_se_month_o
foreach i of numlist 1/22 {
	local j=`i'+3
	gen o_coef_`j'=_b_month_o`i'
	gen o_se_`j'=_se_month_o`i'
}
foreach i of numlist 1/25 {
	gen o_ci_lb`i'=o_coef_`i'-1.96*o_se_`i'
	gen o_ci_ub`i'=o_coef_`i'+1.96*o_se_`i'
}
* organize into format for plot
drop _b* _se*
gen i=1
reshape long o_coef_ o_se_ o_ci_lb o_ci_ub, i(i) j(month)
egen month_since_o=seq(), from (-2) to (22)
order month_since_o, after(i)
* omit binned endpoints
drop if month_since_o==-23 | month_since_o==22
// event study graph 
twoway  (rarea o_ci_lb o_ci_ub month_since_o, cmissing(n n) sort color(gs0%20)) ///
		(line o_coef_ month_since_o, lcolor(black)) ///
		, yline(0, lcolor(black)) xline(0, lpattern(dash)) ///
		graphregion(fcolor(white)) ylabel(-30(10)30, nogrid) xlabel(-20(10)20) ///
		xtitle("Months since Opower Assignment") ytitle("kWh per Month") ///
		legend(order(2 1) label(1 "95% Confidence Interval") label(2 "Event Study Coefficient")) ///
		legend(size(small)) name(F3_wave3, replace)
graph export .\output\figures\F3_wave3.pdf, replace


// wave6 graph------------------------------------------------------------------
use if opower_paper_wave==6 using .\intermediate_data\elec_wave367_winterpre.dta, clear

*bundle event month outside of [-22,21] range
replace month_since_opower=-23 if month_since_opower<-22
replace month_since_opower=22 if month_since_opower>21
* generate event time dummy
tab month_since_opower, gen (m_o)
* generate opower_event time interaction dummy
foreach i of numlist 1/44 {
	gen int_month_since_opower`i'=m_o`i'*opower 
}
* rename all event month before opower as month_o_i
foreach i of numlist 1/23 {
	local j=24-`i'
	rename int_month_since_opower`i' month_o_`j' 
}
* rename all event month after opower as month_oi
foreach i of numlist 25/44 {
	local j=`i'-24
	rename int_month_since_opower`i' month_o`j'
}
* rename start month of opower as month_o
rename int_month_since_opower24 month_o

// run regression
statsby _b _se, clear: reghdfe tot_kwh month_o_23-month_o_2 month_o month_o1-month_o20 $b $c, ///
	absorb($waveXmonth) vce(cluster ky_ba)

* generate coefficient and confidence intervals for each event month
foreach i of numlist 2/23 {
	local j=24-`i'
	gen o_coef_`j'=_b_month_o_`i'
	gen o_se_`j'=_se_month_o_`i'
}
gen o_coef_24=_b_month_o
gen o_se_24=_se_month_o
gen o_coef_23=0
gen o_se_23=0
foreach i of numlist 1/20 {
	local j=`i'+24
	gen o_coef_`j'=_b_month_o`i'
	gen o_se_`j'=_se_month_o`i'
}
foreach i of numlist 1/44 {
	gen o_ci_lb`i'=o_coef_`i'-1.96*o_se_`i'
	gen o_ci_ub`i'=o_coef_`i'+1.96*o_se_`i'
}
drop _b* _se*
gen i=1
reshape long o_coef_ o_se_ o_ci_lb o_ci_ub, i(i) j(month)
egen month_since_o=seq(), from (-23) to (20)
order month_since_o, after(i)
* omit binned endpoints
drop if month_since_o==-23 | month_since_o==22
// event study graph 
twoway  (rarea o_ci_lb o_ci_ub month_since_o, cmissing(n n) sort color(gs0%20)) ///
		(line o_coef_ month_since_o, lcolor(black)) ///
		, yline(0, lcolor(black)) xline(0, lpattern(dash)) ///
		graphregion(fcolor(white)) ylabel(-30(10)30, nogrid) xlabel(-20(10)20) ///
		xtitle("Months since Opower Assignment") ytitle("kWh per Month") ///
		legend(order(2 1) label(1 "95% Confidence Interval") label(2 "Event Study Coefficient")) ///
		legend(size(small)) name(F3_wave6, replace)
graph export .\output\figures\F3_wave6.pdf, replace


// wave7 graph------------------------------------------------------------------
use if opower_paper_wave==7 using .\intermediate_data\elec_wave367_winterpre.dta, clear

*bundle event month outside of [-22,21] range
replace month_since_opower=-23 if month_since_opower<-22
replace month_since_opower=22 if month_since_opower>21
* generate event time dummy
tab month_since_opower, gen (m_o)
* generate opower_event time interaction dummy
foreach i of numlist 1/37 {
	gen int_month_since_opower`i'=m_o`i'*opower 
}
* rename all event month before opower as month_o_i
foreach i of numlist 1/23 {
	local j=24-`i'
	rename int_month_since_opower`i' month_o_`j' 
}
* rename all event month after opower as month_oi
foreach i of numlist 25/37 {
	local j=`i'-24
	rename int_month_since_opower`i' month_o`j'
}
* rename start month of opower as month_o
rename int_month_since_opower24 month_o

// run regression
statsby _b _se, clear: reghdfe tot_kwh month_o_23-month_o_2 month_o month_o1-month_o13 $b $c, ///
	absorb($waveXmonth) vce(cluster ky_ba)

* generate coefficient and confidence intervals for each event month
foreach i of numlist 2/23 {
	local j=24-`i'
	gen o_coef_`j'=_b_month_o_`i'
	gen o_se_`j'=_se_month_o_`i'
}
gen o_coef_24=_b_month_o
gen o_se_24=_se_month_o
gen o_coef_23=0
gen o_se_23=0
foreach i of numlist 1/13 {
	local j=`i'+24
	gen o_coef_`j'=_b_month_o`i'
	gen o_se_`j'=_se_month_o`i'
}
foreach i of numlist 1/37 {
	gen o_ci_lb`i'=o_coef_`i'-1.96*o_se_`i'
	gen o_ci_ub`i'=o_coef_`i'+1.96*o_se_`i'
}
drop _b* _se*
gen i=1
reshape long o_coef_ o_se_ o_ci_lb o_ci_ub, i(i) j(month)
egen month_since_o=seq(), from (-23) to (13)
order month_since_o, after(i)
* omit binned endpoints
drop if month_since_o==-23 | month_since_o==22
// event study graph 
twoway  (rarea o_ci_lb o_ci_ub month_since_o, cmissing(n n) sort color(gs0%20)) ///
		(line o_coef_ month_since_o, lcolor(black)) ///
		, yline(0, lcolor(black)) xline(0, lpattern(dash)) ///
		graphregion(fcolor(white)) ylabel(-30(10)30, nogrid) xlabel(-20(10)20) ///
		xtitle("Months since Opower Assignment") ytitle("kWh per Month") ///
		legend(order(2 1) label(1 "95% Confidence Interval") label(2 "Event Study Coefficient")) ///
		legend(size(small)) name(F3_wave7, replace)
graph export .\output\figures\F3_wave7.pdf, replace

********************************************************************************
capture graph close
capture log close
exit
